Here are the details of the Bayesian multi-level modelling. Our general approach is:
Define Priors
In this section we will specify some priors. We then then use a prior-predictive check to assess whether the prior is reasonable or not (i.e., on the same order of magnitude as our measurements).
Fixed Effects
Our independent variable is a categorical factor with 16 levels. We will drop the intercept from our model and instead fit a coefficent for each factor level (\(y \sim x - 0\)). As our dependant variable will be log-transformed, we can use the priors below:
prior <- c(
set_prior("normal(0,2)", class = "b"),
set_prior("cauchy(0,2)", class = "sigma"))
Group-level Effects
We will keep the default weakly informative priors for the group-level (‘random’) effects. From the brms documentation:
[…] restricted to be non-negative and, by default, have a half student-t prior with 3 degrees of freedom and a scale parameter that depends on the standard deviation of the response after applying the link function. Minimally, the scale parameter is 10. This prior is used (a) to be only very weakly informative in order to influence results as few as possible, while (b) providing at least some regularization to considerably improve convergence and sampling efficiency.
Prior Predictive Check
Now we can specify our Bayesian multi-level model and priors. Note that as we are using sample_prior = 'only', the model will not learn anything from our data.
m_prior <- brm(data = d_eeg,
rms ~ wg-1 + (1|subject),
family = "lognormal",
prior = prior,
iter = n_iter ,
sample_prior = 'only')
We can use this model to generate data.
## Observations: 320,000
## Variables: 2
## $ key <chr> "P2", "P2", "P2", "P2", "P2", "P2", "P2", "P2", "P2", "P2"…
## $ value <dbl> 2.58726956, 0.39915973, 0.60307342, 0.70885249, 0.11650792…
We can see that i) our priors are relatively weak as the predictions span several orders of magnitide, and ii) our empirical data falls within this range.
Compute Posterior
Fit Model to EEG Data
We will now fit the model to the data.
## Family: lognormal
## Links: mu = identity; sigma = identity
## Formula: rms ~ wg - 1 + (1 | subject)
## Data: d_eeg (Number of observations: 400)
## Samples: 4 chains, each with iter = 10000; warmup = 5000; thin = 1;
## total post-warmup samples = 20000
##
## Group-Level Effects:
## ~subject (Number of levels: 25)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.38 0.06 0.28 0.52 1.00 1417 2504
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## wgCM -0.49 0.09 -0.66 -0.31 1.00 877 1502
## wgCMM -0.13 0.09 -0.29 0.05 1.00 845 1513
## wgP2 -0.95 0.09 -1.12 -0.77 1.00 873 1536
## wgP3 -0.82 0.09 -0.99 -0.64 1.01 854 1517
## wgP31M -0.44 0.09 -0.61 -0.26 1.00 861 1337
## wgP3M1 -0.14 0.09 -0.31 0.04 1.00 870 1559
## wgP4 -0.48 0.09 -0.66 -0.31 1.00 884 1412
## wgP4G -0.26 0.09 -0.43 -0.08 1.00 842 1483
## wgP4M 0.15 0.09 -0.02 0.33 1.00 896 1618
## wgP6 -0.48 0.09 -0.66 -0.31 1.00 847 1536
## wgP6M -0.05 0.09 -0.22 0.13 1.00 860 1591
## wgPG -0.93 0.09 -1.11 -0.75 1.00 849 1471
## wgPGG -0.73 0.09 -0.90 -0.55 1.00 851 1504
## wgPM -0.59 0.09 -0.77 -0.42 1.00 852 1373
## wgPMG -0.26 0.09 -0.43 -0.08 1.00 847 1549
## wgPMM -0.07 0.09 -0.24 0.11 1.01 860 1378
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.23 0.01 0.21 0.25 1.00 8172 10803
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
We will now look at the model’s predicts for the average participant (i.e, ignoring the random intercepts).
Fit Model to Psychophysical Data
We will now fit the model to the data.
## Family: lognormal
## Links: mu = identity; sigma = identity
## Formula: threshold ~ wg - 1 + (1 | subject)
## Data: d_dispthresh (Number of observations: 186)
## Samples: 4 chains, each with iter = 10000; warmup = 5000; thin = 1;
## total post-warmup samples = 20000
##
## Group-Level Effects:
## ~subject (Number of levels: 12)
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sd(Intercept) 0.39 0.11 0.23 0.66 1.00 3568 5476
##
## Population-Level Effects:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## wgCM -0.22 0.17 -0.53 0.12 1.00 2280 3879
## wgCMM -1.15 0.17 -1.47 -0.81 1.00 2259 3668
## wgP2 -0.29 0.17 -0.61 0.06 1.00 2578 4062
## wgP3 -0.69 0.17 -1.00 -0.35 1.00 2289 3760
## wgP31M -1.18 0.17 -1.49 -0.85 1.00 2338 4052
## wgP3M1 -1.35 0.17 -1.66 -1.01 1.00 2273 3518
## wgP4 -0.89 0.17 -1.21 -0.55 1.00 2361 3303
## wgP4G -1.22 0.17 -1.54 -0.88 1.00 2484 4249
## wgP4M -1.29 0.17 -1.61 -0.95 1.00 2181 3222
## wgP6 -1.20 0.17 -1.53 -0.86 1.00 2455 3616
## wgP6M -1.42 0.17 -1.74 -1.07 1.00 2552 3730
## wgPG 0.36 0.17 0.04 0.71 1.00 2486 3831
## wgPGG -0.31 0.16 -0.62 0.03 1.00 2381 3981
## wgPM -0.79 0.17 -1.11 -0.44 1.00 2546 3953
## wgPMG -0.96 0.16 -1.27 -0.63 1.00 2452 3772
## wgPMM -1.17 0.17 -1.49 -0.83 1.00 2311 3837
##
## Family Specific Parameters:
## Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
## sigma 0.42 0.02 0.37 0.46 1.00 15742 14912
##
## Samples were drawn using sampling(NUTS). For each parameter, Bulk_ESS
## and Tail_ESS are effective sample size measures, and Rhat is the potential
## scale reduction factor on split chains (at convergence, Rhat = 1).
EEG Control Data
We will also fit models to the control data. As we can see from Figure 2.4, the group differences are much smaller.
## Warning: Bulk Effective Samples Size (ESS) is too low, indicating posterior means and medians may be unreliable.
## Running the chains for more iterations may help. See
## http://mc-stan.org/misc/warnings.html#bulk-ess